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The propagation of electromagnetic waves in isotropic dielectric media with local dispersion is 
studied under the assumption of small but nonvanishing X/£, where A is the wavelength, and £ is 
the characteristic inhomogeneity scale. It is commonly known that, due to nonzero X/£, such waves 
can experience polarization-driven bending of ray trajectories and polarization dynamics that can 
be interpreted as the precession of the wave “spin”. The present work reports how Lagrangians 
describing these effects can be deduced, rather than guessed, within a strictly classical theory. 

In addition to the commonly known ray Lagrangian that features the Berry connection, a simple 
alternative Lagrangian is proposed that naturally has a canonical form. The presented theory 
captures not only the eigenray dynamics but also the dynamics of continuous wave fields and rays 
with mixed polarization, or “entangled” waves. The calculation assumes stationary lossless media 
with isotropic local dispersion, but generalizations to other media are straightforward to do. 


I. INTRODUCTION 


Electromagnetic (EM) waves propagating in inhomo¬ 
geneous linear media exhibit a variety of intriguing phe¬ 
nomena associated with wave polarization. Those include 
precession of the polarization vector, known as the Ry- 
tov rotation [DSl, and also the polarization-driven bend¬ 
ing of ray trajectories, known as the optical Magnus ef¬ 
fect or the optical Hall effect [5HH]- As overviewed re¬ 
cently in Ref. [5], these phenomena can be attributed 
as manifestations of the Berry phase cn], which is a 
fundamental concept emerging also in many other areas 
of physics [IIHS]. Hence the subject has been attract¬ 
ing an increased amount of attention, especially with fo¬ 
cus on variational formulations, which are particularly 
elucidating at studying EM polarization effects in the 
geometrical-optics (GO) limit [71 [T^HTT] . Yet a straight¬ 
forward universal theory is still lacking. Within existing 
models, ray Lagrangians and Hamiltonians either have 
to be guessed or are derived using additional postulates 
like quantization, which is, by definition, alien to classi¬ 
cal theory. Such ad hoc approaches result in a number of 
limitations; e.g., they are not readily applicable to me¬ 
dia with nonlocal and strongly anisotropic dispersion like 
magnetized plasma. Thus, there is a compelling need for 
a generalized and simplified description of EM polariza¬ 
tion effects from first principles. 

The purpose of this paper is to present such a descrip¬ 
tion. More specifically, what we report here is a new 
application of the general theory that we derived earlier 
in Ref. [THl for the Dirac electron. In the present paper, 
we explain how exactly the same theory can be applied 
to GO EM waves and polarization effects in particular. 
As opposed to Ref. [TS], our derivation is entirely classi¬ 
cal. It is also straightforward and elementary in the sense 
that the wave Lagrangian does not need to be guessed, 
as in Refs. [3 US], but is rather deduced (basically, us¬ 
ing nothing more than matrix multiplication) according 
to a formalized algorithm. The known results are suc¬ 
cessfully reproduced and are extended as follows: (i) We 
show that, in addition to the commonly known ray La¬ 


grangian that features the Berry connection [T^2Ij . an 
alternative ray Lagrangian is possible that naturally has 
a canonical form and corresponds to a Hamiltonian sim¬ 
pler than that proposed in Ref. [7]. We explain how the 
two Lagrangians are related and demonstrate their (ap¬ 
proximate) equivalence numerically, (ii) As opposed to 
the result of Ref. [T5] , our Lagrangians are expressed in 
terms of physical time and wave vector, so they capture 
the complete ray dynamics rather than just the ray tra¬ 
jectory. (iii) In addition to eigenray equations, we also 
derive equations for continuous wave fields and rays with 
mixed polarization, or “entangled” waves. This descrip¬ 
tion captures both the Rytov rotation and the optical 
Hall effect simultaneously, and it is also manifestly con¬ 
servative, since the amplitude equations too are derived 
from a variational principle. 


The calculation presented below applies to arbitrary 
linear stationary lossless media with isotropic local dis¬ 
persion, i.e., media whose dielectric and magnetic per¬ 
mittivities are real scalars depending only on spatial co¬ 
ordinates. However, within our theory, generalizations 
to strongly anisotropic media with nonlocal dispersion, 
such as magnetized plasma, are straightforward to do. 
The present paper is intended as an introduction to such 
calculations, which will be reported separately. Thus, be¬ 
low we intentionally focus on a relatively simple problem 
to show how our general theory fits existing literature. 


This work is organized as follows. In Sec. [TJ we de¬ 
fine the basic notation. In Sec. |HI[ we present a La¬ 
grangian formalism to describe EM wave propagation 
in isotropic stationary dielectric media with local dis¬ 
persion. In Sec. IIYI we obtain a reduced model which 
captures first-order polarization effects in transverse EM 
waves. In Sec. [^ we discuss a fluid Lagrangian model, 
which describes the dynamics of the wave envelope, ray 
trajectory, and polarization. In Sec. jvji we derive ray 
equations in canonical form. In Sec. VH the noncanoni- 
cal ray equations are presented and compared with those 
previously reported. In Sec. |vhti we summarize our 
main results. 
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II. NOTATION 

The following notation is used throughout this work. 
The symbol “=” denotes definitions, “h. c.” denotes 
“Hermitian conjugate”, and “c. c.” denotes “complex 
conjugate”. Also, In denotes a unit N x N matrix, 
and hat (') is reserved for differential operators. The 
Minkowski metric is adopted with signature (+, —, —, —), 
so, in particular, d'^x = dt d^x. Generalizations to curved 
metrics are straightforward to apply [52]. Greek indices 
span from 0 to 3 and refer to spacetime coordinates, x^, 
with x° corresponding to the time variable, t; in par¬ 
ticular, dfj, = djdx^. Latin indices span from 1 to 3 
and denote the spatial variables, x* (except where spec- 
ihed otherwise); in particular, di = djdx^. Summation 
over repeated spatial indices is assumed. Also, in Euler- 
Lagrange equations (ELEs), the denotation “Sa :” means, 
as usual, that the corresponding equation was obtained 
by varying the action integral with respect to a. 


III. BASIC EQUATIONS 

A. Photon wave function 

In this paper, we consider EM wave propagation in 
isotropic lossless dielectric media with local linear dis¬ 
persion. In this case, the governing equations for the 
electric field E(t,x) and magnetic field H(t,x) are 

5iE=^VxH, (1) 

5iH = --VxE, (2) 

where e = e(x) is the electric permittivity, and /x = /i(x) 
is the magnetic permeability. 

Let us introduce the following normalized fields. 


E = v^E, (3) 

H = v7^H, (4) 

so one can rewrite Eqs. Q and ([^ as 

atE=-V X H--V(lnV)I) X H, (5) 

n n 

5(H = --VxE-b-V(lnVe) xE, (6) 

n n 


where n(x) = is the refraction index. Due to the 
linearity of these equations, we may formally extend them 
to complex fields. Then, we can express Eqs. ([^ and ®> 
as a vector Schrodinger equation: 

idtip = (7) 

where 



is a six-component wave function, which can be inter¬ 
preted as the photon wave function |22j [23] . The Hamil¬ 
tonian operator is a 6 x 6 matrix given by 

H{x,k) = -X-k + A, (9) 

n 

where k = —zV is the wave vector (momentum) opera¬ 
tor, the A matrices are 6x6 Hermitian matrices, 

^ ^ (-L i) • (“) 


and the matrix A(x) is defined as follows: 


M(x) = - 


n I a • V (In y/e) 


-OL ■ V(lny^^) 


0 


( 11 ) 


Here, the a matrices are 3x3 traceless Hermitian ma¬ 
trices given by [5^ 


/O 0 0\ 

0 0 -z , 

\0 z 0 / 

0 0 z\ 

0 0 0 , 

-z 0 0/ 



a 


Z 



( 12 ) 


(13) 


(14) 


These matrices serve as generators for the vector product. 
Namely, for any two column vectors A and B, one has 

(a • A)B = iA X B, (15) 

A^a^B = -i(A X B)^ (16) 


where T denotes the matrix transpose. 

Let us also point out that, since the Hamiltonian is 
Hermitian, it conserves the wave action: = 0, 

where = (i9t, V), and = (l6,cA/n). More explic¬ 
itly, this can be written as 


d_ / |Ep + |Hp \ ^ rc Re(r xH) - 

dt \ IOtt J n Stt 


(17) 


where the Stt factor was introduced to emphasize the con¬ 
nection to the well known Poynting theorem |25] . 


B. Lagrangian density 


Equation Q has a form akin to the so-called multi- 
symplectic form [26H28] and can be readily given a vari¬ 
ational interpretation. Specifically, Eq. Q can be ob¬ 
tained as the Euler-Lagrange equation, 6A = 0, for the 
action integral 


A = 


-0(t,x) = 


£d^x. 


(18) 
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where the Lagrangian density £ is given by 

£ = 1 (19) 

Here we adopted natural units, such that c = 1. We also 
introduced 

= (""I 

where Z{'x.) is the impedance of the medium, 


that correspond to a local wave vector k. At vanishingly 
small e, these eigenfrequencies are found from the GO 
limit of Eq. namely, 

= i?o(x,k)^, (26) 


where Hq{x, k) = A • k/n. Notice that the matrix A does 
not enter Hq, because A is of the first order in e. 

In Eq. (26), there exists in general six independent 
eigenmodes, since Hq is Hermitian. These modes can be 
readily obtained from 


ZiK) = 



( 21 ) 




1 /(a • k)^ 
0 


0 ^ 

(«-k)V 




(27) 


It is to be noted that no approximations have been used 
in order to obtain this Lagrangian model. Notice also 
that the Lagrangian density (IW ) involves only first-order 
differential operators. In Ref. |i 8 |. we showed that such 
a simple structure is convenient for studying a point- 
particle model of the spin-orbit coupling for the Dirac 
electron. Below, we report how the calculation can also 
be extended to study the effect of polarization of classical 
waves in dielectrics of the specified type. 


Denoting the components of ^ as = (af a^), the 
equation for the first three components is 

ui^ai =^(a ■ k)(Q: • k)ai 
=-^k X (k X ai) 

= ^ [fc^ai - (k • ai)k] , (28) 


IV. REDUCED MODEL 

We consider waves such that 


= ( 22 ) 

where 0{t,x) is some rapid real phase, and ^(t, x) is a 
vector evolving slowly compared to 9(t,x). Hence, we 
assume that the inhomogeneity scale £ is large compared 
to the wavelength; i.e., 

(23) 

where k = |k|, and 

k = Ve (24) 

is the wave vector. Below, we construct a reduced model 
of EM wave propagation using the smallness of e. Unlike 
the standard geometrical-optics (GO) theory, which cor¬ 
responds to a Lagrangian accurate to the zeroth order in 
e HSIISD], this reduced model will yield a ray Lagrangian 
accurate to the first order in e. 


A. Eigenmodes in the limit of vanishing e 

In general, there exist multiple eigenfrequencies, 

w = -dt0, (25) 


where, in the second line, we have used Eq. (15). As 
usual, two eigenmodes are related to the propagation of 
longitudinal (k x ai = 0 ) modes with zero frequency, 
while the other four eigenmodes correspond to the prop¬ 
agation of transversal (k • ai = 0) EM wave modes in the 
limit of vanishing e. 

We will particularly be interested in transversal EM 
modes with positive phase velocities, such that 


uj = k/n. 


(29) 


Corresponding to this frequency, there are two orthonor¬ 
mal eigenvectors hq., which are given by 


where ei(k) and 62 (k) are any two orthonormal vectors 
in the plane normal to et = k/fc. A right-hand conven¬ 
tion is adopted, such that ei x 02 = e^. Obviously, hi 
and ^2 are eigenvectors of Hq. Eor example, 


Uohi — — (k • \)hi 
n 

^ i f (a • k)e2 \ 

nV2 \-{a-k)ei) 

_ k fei 
nV2 V®2 

=ujhi, (31) 


where in the third line, Eq. (151 was used, 
calculation follows for /i 2 . 


A similar 
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B. Eigenmode decomposition 


C. Stern-Gerlach Hamiltonian 


Since hq form a complete basis, one can write ^ = hq(j)‘^, 
where (j)‘^ are scalar functions. We will assume that only 
transverse modes with positive frequencies are actually 
excited (we call these modes “active”), whereas others 
acquire nonzero amplitudes only through the medium in¬ 
homogeneity or the finite width of the EM wave packet 
(we call these modes “passive”). In this case, 


r O(e 0 ), g=l, 2 , 

\ 0{^), q^l,2. 


(32) 


[For future references, note that 0^’^(t,x) describe the 
envelopes corresponding to the linearly polarized modes 
with the electric field E aligned parallel to the unit vec¬ 
tors ei and 62 , respectively.] For given one can 

always calculate the amplitudes of passive modes using 
the complete set of Maxwell’s equations, but this will not 
be necessary for our purposes. As shown in Ref. [18], due 
to the mutual orthogonality of all hq, the contribution of 
passive modes to £ is o(e), so it can be neglected entirely. 
In other words, for the purpose of calculating £, it is suf¬ 
ficient to adopt ^ ~ hi(l)^ + h 24 P‘, even though the true ^ 
may have nonzero projections also on other hq. 

It is convenient to rewrite this decomposition in a ma- 


trix form. 




^ = E(j), 

(33) 

where 




= (t2^ . 

(34) 


Let us search for a tractable expression for U, which 
in Refs. [HI EI] was called the “Stern-Gerlach” Hamilto¬ 
nian. Regarding the first term on the right hand side of 
Eq. (39), a straightforward calculation gives 


S'I'AIS 


_ 1 / 0 a • V In Z\ fei 62 \ 

~2V2n^ l^a-VlnZ 0 j 1^02 -ej 

_ i /e^ \ /VZ x 62 —VZ x ei\ 

~AnZ -e^y l^VZ x ei VZ x 82 j 

= 0. (41) 


Here, we introduced the notation = ef and = e|’. 
Notably, these can be understood as adjoint basis vectors; 
e.g., e^ei = (e^)i(ei)® = 82-81 =0. Hence, Eq. (39) is 
simplified down to 


Z^ = ^[sV(9i.S)-h.c.]. (42) 

Using that = {d'E./dkj){d^kj), one can rewrite U as 
follows [18]: 

Z^=^k-[st(5kS)-(akSt)S]. (43) 

Here, k(t) is given by the lowest-order (in e) ray equation, 


k(t) = -5xw(x,k), 


(44) 


and uj is given by Eq. (29). Using Eq. (35), we obtain 


and 2 is a 6 X 2 matrix having hq as its columns: i.e.. 


:(k) = 




^2 V®2 -ei J ' 


(35) 


[Below, we also consider 2 as a function of in the 
sense that 2 = 2(k(a:^)).] Then, inserting Eq. (33) into 
Eq. ( [I^ , one obtains [IS] 

Z = K.-(t)^ {S-U)(j) + o{e), (36) 


where 


^ ^ Kf/*) - C.C.] , (37) 

E = dt6+^, (38) 

n 

U = 2tM2 + ^ [St7'^(9^2) - c. c.] . (39) 


Here, dt = 9t-l-Vo-V is a convective derivative associated 
to the zeroth order (in e) geometrical optics velocity field, 


. doj 1 
Vo = = -ek. 

ok n 


(40) 


The terms JC and U, which are of order e, represent cor¬ 
rections to the standard, lowest-order GO Lagrangian 
density. Below, U will be calculated explicitly, whereas 
the higher-order terms, o(e), will be neglected. 




dE _ 

( 


d 






1 d 2 9 

0 


Since 8^82 = 0, Eq. (43) leads to 

U = —k • A(k)cry, 
where ay is the y-component of the Pauli matrices, 


(45) 


(Ty - 


0 -i 
i 0 


and A(k) is a vector with components given by 

For example, one may choose [3S] 

kx kz 

kukz 


/ kx kz \ 


81 (k) = 


y/kl+kl J 


, e2(k) = 


/_ ky \ 

kx 

y/W^y 

V 0 y 


(46) 


(47) 


(48) 
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so that 


V. CONTINUOUS WAVE MODEL 


A(k) = 


e_L X Ok 


(49) 


where = (fca;,fcy,0), and ej^ = kj^/fc_L; or, more ex¬ 
plicitly, 


k I 

A(k) = j^\-k: 


(50) 


-L V 0 


D. Lagrangian density: summary 


Lagrangian density as 


£ = -0M dt9 + - 


Substituting Eq. (45) into Eq. (36), we can express the 


+ c.c.) - (j)'^ay'E(j), (51) 


where 

E = k-A(k). 

Let us also use a variable transformation, 
())(t,x) = Qr]{t,x), 

where 

-0' 


(52) 

(53) 

(54) 


and ? 7 (<,x) is a new vector with components denoted as 
^(^:X)=(^+). (55) 


Hence, the Lagrangian density (51) can be expressed as 

£ = —T]^ {dtO + k/n) T] 

+ ^{il^dtr] - c.c.) - r]^a^T,r], (56) 

where is another Pauli matrix, 

'1 0 


a, = 


0 -1 


(57) 


Here, ri± (t, x) describe envelopes corresponding to right- 
hand and left-hand circularly polarized modes, respec¬ 
tively (as defined from the point of view of the source). 
The first line of the right hand side of Eq. (56) repre¬ 


sents the lowest order GO Lagrangian density for these 
modes. Additionally, the second line of Eq. (56), which 


contains 0{e) terms, introduces polarization effects, as 
will be explained below. 


Substituting Eqs. (29) and (44) in Eq. (52), we can 
rewrite £ as 


£ = -uM dtO + - 


+ ^iv^dt-T] - c.c.) - r]'^az^{x,k)r], (58) 


where 


E(x,k) = 4jA(k) • Vn. 


(59) 


In this form, the Lagrangian density is analogous to that 
of a semiclassical Pauli particle and thus can be ap¬ 
proached similarly [3T]. Let us adopt the representation 
?7 = z\/X^ where I(t,x) = is the action density, and 
z{t,x) is a unit polarization vector {z^z = 1). Since the 
common phase of the two components of z can be at¬ 
tributed to 6, we can parameterize z in terms of just two 
real functions, 'd(t,x) and 

e*’’ cos(C/2) \ 

\^e-*’’sin(C/2)j ■ 

Like in the case of the Pauli particle, C determines the 
relative fraction of “spin-up” and “spin-down” quanta, 
i.e., those corresponding to left-hand and right-hand po¬ 
larizations. Also, one can understand S = zf(Tz/2 as 
the wave “spin vector” (cr denotes the three Pauli matri¬ 
ces, as usual [Sin]), or, up to a constant factor, as the 
Stokes vector [321131] ■ 

Using the above parameterization, the Lagrangian den¬ 
sity is rewritten as 


£ = -I 


dtO -\ -1- {dtd + E) cos C 

n 


(60) 


which leads to four ELEs. The first one is the action 
conservation theorem, 

66: dtl + V ■ (IV) = 0. (61) 

The flow velocity is given by V = Vq -I- u, and 


d /k-Vd 


dk 


nk 


+ E I cos C- 


(62) 


Notice that u represents the polarization-driven deflec¬ 
tion of the ray’s “center of gravity” predicted in Refs. [3- 
[7]. The second ELE is a Hamilton-Jacobi equation, 

61: dtO + — + {dtd + E) cos C = 0, (63) 

n 

whose gradient yields an equation for dtk, i.e., the mo¬ 
mentum equation [31j . The third ELE is 


: dtd = —E. 


(64) 
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As it will become clear below, this describes the rotation 
of the wave polarization. Finally, the fourth ELE is 


Si9 : 9((XcosC) + V • (Ivq cosC) = 0. (65) 


Together, Eqs. (61)-(65) provide a complete “fluid” de¬ 


scription of continuous waves. Note that waves are al¬ 
lowed to be nonstationary and “entangled”, i.e., contain 
mixed polarization. In fact, a combination of Eq. (65) 


with Eq. (61) gives 


9((cosC) + V • V(cosC) = • (XucosC), ( 66 ) 


which shows that ( is generally not conserved; i.e., the 
numbers of “spin-up” and “spin-down” quanta necessar¬ 
ily oscillate, unless a wave is homogeneous. This can be 
interpreted as an effective zitterbewegung of a classical 
wave in an inhomogeneous medium. 


B. Polarization dynamics 


To better understand the polarization equations, let 
us rewrite Eq. (0 as an equation in the basis of linearly 
polarized modes, i.e., for $(t) = 

4) = QZ = —iT,QazZ = = —iT,ay^. 

(73) 


[This equation could also be obtained if the ray equations 
were derived directly from the Lagrangian density (51).] 
Since E is a scalar, and ay is constant, this can be readily 
integrated, yielding [3S] 

$ = exp(—i0tTj,)<i)o = (I 2 cos 0 — iay sin 0)$o, (74) 


where 0 (t) = S(F) df can be recognized as the wave 
Berry phase |TU], E(t) = E(X(t), P(t)) and $0 is a vector 
determined by initial conditions. This result can be also 
be expressed explicitly as follows: 


VI. RAY DYNAMICS: CANONICAL 
REPRESENTATION 


$ = 


cos 0 — sin 0 \ 
sin 0 cos 0 y ° 


(75) 


A. Basic equations 


The ray equations corresponding to the above field 
equations can be obtained as a point-particle limit. In 
this limit, I can be approximated with a delta function 

X(t,x) =(5(x-X(t)), (67) 


where X is the location of the wave packet. As in 


Refs. [HISI], the Lagrangian density ( [5^ yields a point- 
particle Lagrangian, L = J £ d^x, specifically, 

L = P -X- — + - (Z^Z - Z^Z) - E(X, P) Z^azZ. 
n 2 

( 68 ) 


It is seen then that the polarization of the EM field ro¬ 
tates at the rate E(t) in the reference frame defined by the 
basis vectors ( 01 , 62 ). Such rotation of the polarization 
plane, also known as Rytov rotation, was studied theo¬ 
retically in Refs. Ulira and observed experimentally in 
Ref. [3]. Clearly, Eq. (64) describes the same effect. 


C. Ray dynamics for pure states 


If a ray corresponds to a strictly circular polarization, 
such that azZ = ±z, the Lagrangian ( 68 ) can be simpli¬ 
fied down to 


L± = P-X- —tS(X,P), 


(76) 


Here P(f) = V0(t,X(t)) is the canonical momentum, 
P = jPj, and Z{t) = z{t,X{t)) is a complex scalar func¬ 
tion. The refraction index and E are evaluated at x = X 
and k = P; e.g., n = n{X{t)). Also, the speed of light 
constant, c, has been reintroduced for clarity. 

Treating X(t), P(t), Z{t), and Z^{t) as independent 
variables, we obtain the following ELEs: 


SP : 

X(t) = ^ + (5pE) Z^azZ, 
nF 

(69) 

SX : 

P(t) = ^Vn-(9xE)ZV,Z, 

(70) 

6Z^ : 

Z(t) = —iZ,azZ, 

(71) 

SZ : 

Z\t) = i£Z^az. 

(72) 


where L± governs the propagation of right and left po¬ 
larization modes, respectively. This Lagrangian has a 
canonical form, L± = P • X — H±{X, P), where H± is a 
Hamiltonian given by 

The variables X and P serve as the canonical coordi¬ 
nate and momentum, so they satisfy canonical Hamil¬ 
ton’s equations, 

SP : X{t) = ^± 5pE, (78) 

nF 

6X: P(t) = ^VnT^xE. (79) 


Together with Eqs. ( |50[ ) and (59), Eqs. (69)-(72) form a 
complete set of equations. The first terms on the right 


hand side of Eqs. (69) and (70) describe the ray dynamics 
in the GO limit. The second terms describe the coupling 
of the mode polarization and the ray curvature. 


Since H± is time-independent, one also readily obtains 
energy (frequency) conservation, 

i7±(X(t),P(t)) = const, (80) 

along the ray trajectory. 
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FIG . 1: Comparison between ray trajectories predicted by the canonical Lagrangian (761 and the non-canonic al L agran gian 
(811. The red line represents X(t) — A(P(t)), where X(t) and P(t) are found by numerical integration of Eqs. (78l and (79l. 
The blue line shows x(t), as found by numerical integration of Eqs. (901 and (911. The black dashed line plots the “spinless” 
ray trajectory governed by the lowest-order GO ray Lagrangian Lq given by Eq. ( |93[ ). The refraction index is chosen to be 
n(x) = 1 + exp[—+ z'^)/l^\, the characteristic width is ^ = 10, and A is chosen in the form (501. The length unit is 
a = ^/2/po- The initial location is xq = (0,10,0), and the initial momentum is Pq = (1,0, l)/a, so e ~ l/(po^) ~ 0.07. It 
is seen that the canonical Lagrangian and the noncanonical Lagrangian predict results that, for the specified parameters, are 
essentially indistinguishable from each other yet differ noticeably from those yielded by Lq. 


VII. RAY DYNAMICS: NON-CANONICAL 
REPRESENTATION 

A. Ray variables 

If the point-particle limit is taken without explicitly 
invoking Eq. (52), an alternative representation of the 


ray Lagrangian can be obtained that connects our results 
to those found in existing literature. For pure states, this 
alternative Lagrangian is given by 

cp 


L± = p • X - 


i(x) 


TP • A(p), 


(81) 


where p • A(p) is known as the “Berry connection” term 
[19]. Notice, in particular, that adding i9kx(k) to A(k), 
where x(k) is an arbitrary scalar function, changes L± 
merely by a perfect time derivative and thus does not 
affect the motion equations. Also notice that we have in¬ 
troduced a different notation for the ray variables [(x, p) 
instead of (X, P)] for the following reason. 

On one hand, as an approximation of the Lagrangian 


(76), Eq. (81) is expected to yield dynamics similar to 


that yielded by Eq. (76). On the other hand, notice that 


Eq. (81) does not have a canonical form. This means 
that (x(p) do not satisfy Hamilton’s equations and thus, 
clearly, cannot be the same as (X,P). To understand 
the connection between the two sets of variables, let us 
rewrite Eq. (81) as 


^ . CD d , . , , dA 

L± = p • X-— =F (P • A ±p • —. (82 

n dt dt 

Dropping the perfect time derivative and introducing 

q = x±A, (83) 


we obtain a Lagrangian in a canonical form: 

cp 


L± = p • q - 


n(q=F A(p))' 


(84) 


The quantity cp/n{'x) serves as the canonical energy and 
is conserved, 


cp/u(x) = const. 


(85) 


Also, since l/u(x) is assumed smooth, one can replace 
this function with its first-order Taylor expansion. Then, 

L±=p-q -^ =F ^^A(p) • Vn(q), (86) 

n(q) n2(q) 

where we omitted terms ©(e^), as usual. By comparing 
this with Eq. (76), we find that 


x = XtA(P), 


p = P. 


(87) 


Notice that |x — X| is of the order of the wavelength, 
i.e., small enough to make x and X equally physical as 
measures of the ray location. 


B. Ray equations 

The equations for x and p can be obtained by com¬ 
bining Eqs. (78), (79), and (87), or they can be derived 


directly as ELEs corresponding to the Lagrangian (81): 


(5p : 

x(t) = — ± p X (Vp X A), 

(88) 


np 


5x : 

pW = 

(89) 






























Using Eq. (50), we can also rewrite them as 


cp P X p 

x(t) = ± 

(90) 

np 


II 

<j 

s 

(91) 


Equations (90) and (91) were reported previously in 
Refs. [HITB]. The equations presented in Refs. IMIITKIIT^ 


m can also be obtained from Eqs. ( [90| and ([9R for the 
rescaled momentum p = np/{cp) and rescaled “time” 
ds = cdt/n{x{t)). Specifically, one gets 


x'=^±A^, 
P P'^ 


p' = Vn, 


(92) 


where primes denote derivatives with respect to s, and 
we used that A = n(x)/p is a constant of motion, as seen 
from Eq. (85 I. 


non-canonical (81) representations. For completeness, we 


also show the calculated ray trajectory as determined by 
the lowest-order GO ray Lagrangian, 


Lo = P • X - cp/n(x). 


(93) 


which does not account for polarization effects. As antic¬ 
ipated, the ray trajectories predicted by the Lagrangians 
(76) and (81) are almost identical and yet differ no¬ 


ticeably from the “spinless” ray trajectory predicted by 
Eq. (93), namely. 


Sp : 

II 

'TS 1^ 

(94) 

5x : 

p(t) = -^Vn. 

(95) 


VIII. CONCLUSIONS 


C. Comparison of the two models 


As we showed explicitly in Sec. VII A| (a nd also by con¬ 
struction), the canonical Lagrangian and the non- 
canonical Lagrangian ( |8l| ) differ only by O(e^) and thus 
are equivalent within their accuracy domain. (The same 
applies to the earlier theories [SlIIlHinillMI] too, where 
the ray Lagrangians are also derived to the first order in 
e.) This means, in particular, that the effect of the Berry 
connection term in Eq. (811 can be attributed simply to 


the choice of coordinates, while the underlying physics 
can also be described by the canonical Lagrangian ([77]). 
The advantage of the Lagrangian (76) is its manifestly 


symplectic structure, which is convenient, for instance, 
for numerical simulations [251. On the other hand, the 


Lagrangian (81) leads to “gauge-invariant” equations, in 


the sense that they are manifestly independent of the 
coordinate representation for the basis vectors, ei and 
02 . Also, assuming that A is chosen in the form (50), 


using the non-canonical Lagrangian (81) can be advan¬ 
tageous when + pI —)• 0, because the right-hand side 
in Eqs. ( [90| and ( |9l| remains finite in this case unless 
p —>■ 0. Therefore, whether the canonical or noncanoni- 
cal form is more convenient in a given case depends on 
the specific application. 

To illustrate how accurate the agreement is between 
the two models, we also performed comparative numer¬ 
ical simulations. Figure 1 shows the ray trajectories for 
a right polarized wave using the canonical (76) and the 


In this paper, we study the propagation of electromag¬ 
netic waves in isotropic dielectric media with local dis¬ 
persion under the assumption of small but nonvanishing 
e = X/i, where A is the wavelength, and £ is the charac¬ 
teristic inhomogeneity scale. It is commonly known that, 
due to nonzero e, such waves can experience polarization- 
driven bending of ray trajectories and polarization dy¬ 
namics that can be interpreted as the precession of the 
wave “spin”. Here, we report how Lagrangians describing 
these effects can be deduced, rather than guessed, within 
a strictly classical theory. In addition to the commonly 
known ray Lagrangian that features the Berry connec¬ 
tion, a simple alternative Lagrangian is also proposed 
that naturally has a canonical form. We explain how 
the two Lagrangians are related and demonstrate their 
equivalence numerically. The presented theory captures 
not only the eigenray dynamics but also the dynamics of 
continuous wave fields and rays with mixed polarization, 
or “entangled” waves. Our calculation assumes station¬ 
ary media with isotropic local dispersion, but generaliza¬ 
tions to other media are straightforward to do, as will be 
reported separately. 
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